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Abstract 

Semiclassical calculations using the Herman-Kluk initial value treatment are performed to deter- 
mine energy eigenvalues of bound and resonance states of the collinear helium atom. Both the eZe 
configuration (where the classical motion is fully chaotic) and the Zee configuration (where the 
classical dynamics is nearly integrable) are treated. The classical motion is regularized to remove 
singularities that occur when the electrons collide with the nucleus. Very good agreement is ob- 
tained with quantum energies for bound and resonance states calculated by the complex rotation 
method. 
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I. INTRODUCTION 



The inability of the old quantum theory to determine the energy spectrum of helium was 
one of the factors motivating the development of modern quantum mechanics. With the 
advent of the modern theory, attempts to treat helium within a classical framework continued 
with the application of semiclassical methods, obtained from an asymptotic analysis of the 
Schrodinger equation, [lj]. These studies eventually identified two key reasons for the failure 
of old theory's attempts to describe the energy levels of this atom: the absence of the Maslov 
index in the older work, [2]] and the largely chaotic nature of the classical electronic dynamics 
in this system. 3), Q] Due to this last factor, truly satisfactory semiclassical treatments had 
to await development of techniques for the quantization of chaotic motion. The first modern 
semiclassical calculations to properly account for the chaotic dynamics in helium appeared 
in the early 199O's.0, Q] These and subsequent works have enabled the quantization of a 
large number of bound and resonance states of the atom. [l|, [(], Q, [s, 9] 

Although other restricted-configuration coplanar studies have been carried out,[lQ!] most 



recent semiclassical treatments are based on a collinear model for the classical motion of t 
electrons, with corrections added for stable motion in the remaining degrees of freedom. 
Collinear helium can exist in two configurations: 



lie 



The eZe configuration, characterized by the two electrons on different sides of the nu- 
cleus. This configuration is relevant for quantization of bound and resonant states of 
the three-dimensional atom with maximal positive values of the Stark quantum num- 
ber K including, apparently, the ground state of the atom. The classical motion for 
this configuration is believed to be completely chaotic, i.e., the periodic orbits are all 
linearly unstable and proliferate exponentially with orbit length. Semiclassical quanti- 



.0,0,9 



zation in this case was performed by a periodic orbit cycle expansion technique 

The Zee configuration, characterized by the two electrons on the same side of the nu- 
cleus. This configuration is implicated in the quantization of certain resonance states 
with maximal negative values of K. In strong contrast to the eZe case, the classical 
behavior for the Zee configuration is nearly integrable, with invariant tori centered 



about a stable periodic orbit describing the so-called "frozen planet" motion. 



121 ] Quan- 



tization for this case was achieved by applying the Einstein-Brillouin-Keller (EBK) 
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treatment. 




The justification for treating helium with reduced dimensionality is that many highly 
excited quantum states are dominated by the classical collinear motion. This conclusion 
rests on three arguments: [l| (a) A scaling analysis shows that the classical dynamics of 
helium with total orbital angular momentum L at energy E is equivalent to the dynamics of 
the system with angular momentum |i?| 1,/2 L at energy —1. Thus, as the energy approaches 
the double ionization limit E = 0, the classical motion relevant for states with low and 
moderate angular momentum quantum numbers becomes that for L = 0, which restricts 
the three particles to a plane, (b) It is known that the collinear eZe and Zee configurations 
form invariant subspaces in the full-dimensional phase space, (c) It can be shown that 
the collinear motion is stable with respect to bending deformations perpendicular to the 
linear configuration. On the basis of these arguments many states of the three-dimensional 
atom are expected to be structured about the backbone of collinear classical motion. This 
reduction in dimensionality was crucial for the ability to quantize states associated with the 
eZe configuration since it enabled a thorough study and coding of the system's periodic 
orbits. 

Although the existing semiclassical treatments quite successfully reproduce energy levels 
for many states of the three-dimensional atom, full-dimensional calculations of helium are 
highly desirable. One problem with the collinear treatments is that they are applicable 
only to states with near-extreme values of the Stark quantum number; they are not truly 
appropriate for more general states. A second limitation of these approaches is they are 
not strictly justified for low- lying excited states of the atom where the arguments favoring 
collinearity, or even coplanarity, break down. This comment applies even to 5* states which 
should, arguably, be treated semiclassically with nonzero L due to the Langer correction. [y| 
One consequence is that the accuracy of the semiclassical approximation for the ground 
state of helium has not yet been definitively established. Unfortunately, the semiclassical 
methods that have been used to quantize the chaotic classical motion for collinear eZe helium 
cannot be easily applied to higher- dimensional systems where the dynamics generally have 
a mixed chaotic-regular nature and where a systematic search and characterization of all 
periodic orbits becomes unfeasible. Clearly, this limitation also blocks extension of these 
semiclassical techniques to atomic and molecular systems having more than two electrons. 
Thus, it is important to develop other semiclassical approaches that are capable of treating 
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higher- dimensional many-electron systems. 

Semiclassical initial value representation (IVR) methods Q, 12, Q, 19] are natural 
candidates for this role since they are capable of semiclassically propagating wave functions 
for systems with many degrees of freedom. Here we test the applicability of one such method 
for the quantization of many-body Coulombic systems by applying it to the collinear helium 
atom. In contrast to much of the previous semiclassical work for this system, no corrections 
for the motion perpendicular to the linear configuration are incorporated in the present 
treatment, so that the results apply to the collinear system only, not to the three-dimensional 
atom. 

Alt 
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lough IVR treatments have been previously applied to single-electron systems |1 8 
the present application to a system with two electrons is nontrivial. The most 
serious difficulty concerns the strongly chaotic nature of the classical motion for the eZe 
configuration. It is well established that such behavior creates numerical convergence diffi 
culties in the propagation of systems by IVR methods beyond very short times. 23l . 
Thus, calculations for times long enough to extract quantization information are expected 
to be challenging. A related problem is associated with the classical instability of the helium 
atom with respect to ionization. In contrast to the quantum system, the classical atom can 
autoionize at all energies, even those below the first ionization limit. Thus, in general, one 
electron of the atom escapes to infinity after a few revolutions, leaving a bound He + ion. 
This is an extreme example of the zero-point energy problem 26|, |27| that plagues classical 
trajectory simulations of quantum evolution. More specifically, the potential difficulty in the 
present context arises from the reliance of IVR methods on general classical trajectories that 
are almost never true periodic orbits. Quantization is obtained from those trajectories that 
shadow the periodic orbits for long time intervals. The short lifetime of many trajectories, 
resulting from the classical autoionization, may require a drastic increase in the total number 
of trajectories that must be calculated. A final problem anticipated for IVR treatments of 
the collinear atom in both configurations (and one that is common to all numerical semiclas- 
sical calculations for these systems) arises because the electrons undergo repeated collisions 
with the nucleus, where the Coulombic interaction becomes infinite. At such collisions, the 
classical equations of motion become singular, and they cannot be integrated by ordinary 
methods. Numerical solution of the classical equations of motion for collinear helium thus 
requires regularization techniques. 
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The remainder of this paper is organized as follows. Section II presents the expressions 
that form the basis of the present treatment. Section III describes the semiclassical calcula- 
tions, Sec. IV presents the numerical results, and Sec. V discusses and summarizes the main 
points of our work. An Appendix presents details of supporting quantum calculations for 
the energy levels of the collinear atom. 



II. THEORY 

In atomic units, the Hamiltonian for the collinear helium atom is given by 

H = A + rA_2__i + i 

2 2 q 1 q 2 |g 2 ±9i| 

where q T = (91,92) are the coordinates and p T = (pi,P2) are the momenta of the two 
electrons. The + and — signs in the last term apply to the eZe and Zee configurations, 
respectively. 

Here we wish to determine the energies of bound and resonance states for these two 
configurations. In a purely quantum calculation, such energies can be extracted from peaks 
in the spectrum 



exp(iEt) c(t) dt 







(2) 



F(E)- 

where c(t) is the autocorrelation function 

c(t) = J dr' J drV* (r')K(r',r,t)%(r). (3) 

Here \l/o denotes an initial state which projects onto energy eigenstates of interest and 

K(r', r, t) = (r'| exp(-iHt) |r) (4) 

is the propagator, describing the time evolution operator in the representation of positions 
r T = (ri,r 2 ). In the present semiclassical calculations, this propagator is approximated 
semiclassically. 

The specific semiclassical treatment used here is that of Herman and Kluk (HK) jlol . 



who expressed the propagator in the IVR form 

i 2 gd : 

(27T) 



K HK (v',r,t)= [[ ^^(r'|q t ,p t ) J R qpt exp(^ qp ,)(q,p|r), (5) 
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where the integrals are over momenta p and coordinates q that serve as initial conditions for 
classical trajectories. The trajectories carry these phase space variables at time to values 
(pt, q t ) at time t. The quantities 



(i'|q-p) = ( ) <^p 



1/4 

(6) 



■-(r - q) T 7(r - q) + ip T '(r - q) 



are Gaussian coherent state functions. 7 is a 2 x 2 diagonal matrix with real nonzero 
elements 71,72 which determine the Gaussian widths along 77 and r 2 . The function S^pt = 
f d£'(p t /q t / — H) is the action along the trajectory and the prefactor i? qP i has the form 



-Rqpt = < det 



1 / dqt , -idpt .dq.t . _i<9p 4 
1 ~5 — 7 — — 7 + ^7 



1/2 

(7) 





(8) 



2 \ <9q <9p <9p <9q 

where each term under the radical is a 2 x 2 block of the 4-dimensional monodromy matrix 
defined by 

<9q f /<9q <9q 4 /<9p 
<9p t /<9q dp t /dp 

A potential difficulty with applying the HK propagator in the present case is that (j3J) 
is technically valid only for systems with wave functions that obey boundary conditions 
at Ti = ±00. In contrast, the wave functions of collinear helium are required to obey 
boundary conditions at r« = (where they must vanish [4]). Due to the Gaussian form 
of the functions (r'|q 4 ,p 4 ), the wave functions calculated using the HK treatment will not 
satisfy this property. To remedy this problem, the IVR expression for the propagator should 
be modified to replace the Gaussian coherent state functions with certain functions that 
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Fortunately, however, 



vanish at Ti = and obey some additional requirements, 
these modifications are not needed for the present calculations. The states \l/o( r ); which 
we use to prepare the system at time and to probe the positions of the particles at 
times t, are localized in regions far from r» = 0. Thus, the initial states effectively obey 
the proper boundary conditions and the function c(t) is insensitive to the behavior of the 
propagated wave functions ^t(r') = J drK HK (r', r, t)\E r ( r ) near r[ = 0. The results of an 
HK calculation should then be the same as for a modified IVR propagator in which the 
coherent state functions are replaced by new functions that are strongly non-Gaussian at 
T\ = (as needed to obey the boundary conditions) but become almost Gaussian far from 
the boundaries. This was confirmed by calculations in which the (r|q, p) and (r'(q t , p t ) were 
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replaced by combinations of Gaussian coherent state functions which vanish at rj = 0. The 
results were essentially identical to those obtained with the ordinary HK expression. 

As mentioned above, it is necessary to regularize the classical equations of motion because 
they become singular when electrons collide with the nucleus (qu = or q 2 t = 0) and cannot 
be integrated by applying ordinary numerical techniques. As a prelude to this regularization, 
it is convenient to transform the coordinates and momenta to scaled variables defined by 

ft = |£|ft, Pi=\E\- l ' 2 Pi (i = l,2) (9) 

where E is the energy of the atom. In terms of these variables, the Hamiltonian can be 
expressed as H = \E\H where 

£ = pi + g_2_2 + 1 

2 2 q 1 q 2 |g 2 ±<?i| 

has an energy-independent value. For the present treatment of bound and resonance states 
(for which E < 0) we set H = — 1. The physical variables (q t , p t ) at time t can be calculated 
by solving Hamilton's equations of motion with Hamiltonian H to obtain the scaled variables 
at the time t, defined by 

t=\Ef 2 t, (11) 

and then applying (jUJ). 

To complete the regularization procedure, we apply a transformation similar to that 



proposed by Kustaanheimo and Stiefel 29j, |30j. In our case, this involves defining a new 



Hamiltonian 7i = q~\q2{H +1) which is numerically equal to zero and is free from singularities. 
Introducing the regularized variables (Qi,Pi), defined via the equations 

qi = Ql ft = ^ (< = 1,2), (12) 

H takes the form 

n = iQlP * g — - + ^) + ^ + j^^fi ■ ( 13 ) 

With this Hamiltonian, the classical equations of motion for dQi T /dr and dPi T /dr are solved 
with initial conditions (Qi, Pi) at time to determine values of (Qi T , Pi T ) as a function of a 
new time variable r. This variable is related to the scaled time t by the differential equation 

f T = QlrQlr, (14) 
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which must be integrated along with Hamilton's equations. The physical variables p it , qu at 
time t are then obtained by reversing the regularization and scaling transformations. 

The physical action integral St needed to form the HK integrand can be obtained from 
the regularized variables by applying 

S t = \E\- l ' 2 S T1 (15) 
where the scaled action S T is determined by integrating the equation 

dr ~ 8 +A«ir + «ar; \Ql±Qiy W 

To evaluate the monodromy elements needed for the semiclassical calculation, we begin 
with the relations 



\E\q it = Q 



2 

IT ' 



2QirPit = V\E\ Pir, (17) 

which are obtained by combining and ffT2"j) . Treating (Qi T , Pj T ) as functions of t and the 
initial conditions a = {qi,q2,Pi,P2}, and taking partial derivatives with respect to any one 
of these variables, aj, we obtain 



oaj \ aaj J t oaj 



Q „ p = - pa to*) + m + * m. (18) 



daij \ dctj J t 2 \ dctj J t 4y/\E\ daj 

To determine the partial derivatives of Qi T and Pi T with respect to aj appearing above, 
we treat Qi T and Pj r as functions of the initial values Qi,Pi and of r. Differentiating with 
respect to <x,- and applying the equations of motion then gives 

dQ iT \ dQ iT dQ k dQ iT dP k dH f dr 



daj J . dQk daj dPk daj dP iT \ da 3 , t 
dP tT \ dP lT dQ k dP lT dP k dH (dr \ 



daj J t dQk daj dPk daj dQ ir \daj 7 , 

Finally, to determine (dr/da) t in the above equations, we differentiate t (considered as a 
function of initial conditions and t) with respect to aj and apply ffTTj) and (TT3I) . This yields 

3 p-^d\E\_ __d^dQk dl^dPk n2 n 2 ( dr \ 
2 VlEl ^ t ~Wk^ + dPkda~ + QlTQ2r \da~j)t 
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from which (dr/da) t can be obtained. These equations thus express the physical monodromy 
matrix in terms of the regularized monodromy matrix and the partial derivatives dt/dQk 
and dt/dPk- The various quantities needed to calculate the physical monodromy matrix can 
now be formed by solving a set of first-order ordinary differential equations in the variable 
r, obtained by differentiating Hamilton's equations (for dQi T /dr and dP^/dr) and (|T4l) (for 
dt/dr) with respect to Qk and P^. 

We mention that an alternative procedure for calculating the physical monodromy matrix 
can be formulated in terms of a transformation to so-called local coordinates j^]. However, 
the present approach was adopted here since numerical tests showed that it allows larger 
time steps for the integration of the differential equations. 

Note that the differential equations that are solved here remain well behaved for all 
values of r. The singularities in the physical variables (monodromy elements and momenta), 
occurring at collisions with the nucleus, arise only from (I18j) - (j20j) . describing transformations 
from the regularized variables. Nevertheless, such singularities cause the prefactor _R qpt to 
behave as q~^ for qn — » 0, raising the suspicion that the HK integral diverges for the 
present system. Fortunately, these singularities in i? qpt are cancelled by factors such as 
exp[— (p it — PioY/^i] appearing in the expressions for c(t) [see ff23l) below], which cause the 

— 1/2 

HK integrand to vanish when q it — > since pn oc q it in such cases. 

While the regularization presented above provides a unique continuation of the classical 
equations of motion when one of the coordinates qu or qit becomes zero (a binary collision), 
no unique continuation exists when both qu and qit become zero simultaneously (a triple 
collision). Fortunately, the initial collisions leading to an exact triple collision form a set of 
zero measure in the system's phase space and do not contribute to calculations of the HK 
propagator. 

III. SEMICLASSICAL CALCULATIONS 

The symmetry of the Hamiltonian in ([T]) implies that two distinct kinds of energy eigen- 
states exist: states that are symmetric with respect to electron exchange (which can be 
identified with singlet spin states) and states that are antisymmetric with respect to elec- 
tron exchange (which can be identified with triplet spin states). As in the case of exact 
quantum evolution, one expects that energy spectra for the two kinds of states can be ob- 
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tained separately by semiclassically propagating initial states of definite exchange symmetry. 
Thus, in our calculations, we choose the functions ^ in (j3J) to have the (anti) symmetrized 
form 

*o(r) = N[(ri\qw,pio)(r 2 \q2o,P2o) ± <r 1 1 g 2 o , P20 > 1 9io , Pio> ] , (21) 

where (ri\qjo,Pjo) are one-dimensional Gaussian coherent state functions involving param- 
eters (qjOjPjo) and Af is a normalization constant. We choose the values of jj for these 
functions to be equal to those used to construct the functions (r|q, p) = (ri\qi, p±) (r 2 \q2, P2) 
in the propagator. However, rather than propagate a combination of coherent state func- 
tions, it is more convenient to propagate a single, unsymmetrized, "initial" state and project 



out the desired symmetry component with a properly symmetrized "final" state. 23j Applied 
to the autocorrelation function this treatment yields 
r r q d^p 

\t) = 2J\f 2 II - 2 [(gi 0) Plo|git,Plt)(920 ) P2ok2t,P2t) ± (?20,^20|git,^lt)(gi0,^10|g2t,P2t)] 



(2 



7T 



Rwt exp (iSqpt) (qi,Pi\ qw , Pw ) (<?2 , P2 \ q20 , P20) ■ (22) 



The product of overlap factors on the extreme right in ( 1221) can be expressed in the 
analytical form 

(q, Plqo, Po) = exp [-|(q - q ) T 7(q - qo) - |(P - Po) T 7 _1 (P - Po)] 

xexp[|(p + p ) T (q-q )]. (23) 

In our calculations, we used the Gaussian distribution, formed from the modulus of this 
expression, to sample points in phase space for a Monte Carlo evaluation of the integrals 
over (q, p) in ( !22l) . 

To calculate the integrand in (122jl . the transformations ([9]) and (TI2I) were applied to de- 
termine values for the regularized variables Qi, Pi for each physical point (q, p) sampled. 
Using these as initial conditions, Hamilton's equations for Q ir and Pj r were solved numeri- 
cally, along with (I14p and (I16p . as well as the r-dependent differential equations for elements 
[e.g. dQi T /dQj] of the regularized monodromy matrix and the quantities dt/dQi, di/dR. 



This set of 26 equations was integrated using a variable stepsize Runge-Kutta method 31]. 
The values of r, corresponding to an evenly-spaced sequence of desired times t, were deter- 
mined by a root-solving procedure and, at such times, the transformations were reversed 
to obtain the physical coordinates, momenta, action, and monodromy elements, needed for 
construction the integrand. 
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The final step in the calculation, namely, the evaluation of the Fourier transform in 



was performed by the harmonic inversion procedure. 



32 



33. 



34[ which effectively fits the 



computed autocorrelation function to an expression of the form 



c HK {t) = J2 d k e M-iE k t), 



(24) 



extracting values for the complex parameters E k and d^. This method was preferred over 
ordinary Fourier transform techniques since it was capable of producing spectra of higher 
resolution spectra from our autocorrelation data, which were calculated for relatively short 
time durations. 

Numerical difficulties, caused by the highly chaotic nature of the classical dynamics for 
the eZe configuration, made it necessary to apply some additional approximations. As is 
well known, the chaotic motion causes the modulus of the HK prefactor -R qp t to increase 



approximately exponentially with time. 
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251 ] This systematic increase is distinct from 



singularities in the prefactor that occur when an electron collides with the nucleus. In 
principle, chaos should also cause the phases in the HK integrand to vary rapidly as a 
function of initial condition and the resulting cancellations should compensate for the growth 
of the prefactor, reducing the autocorrelation function to physically correct values. However, 
as time progresses, the prefactor becomes so large that the number of Monte Carlo-sampled 
trajectories needed for this cancellation exceeds practical limits. As a result, the calculation 
cannot be made to converge beyond a short time. The practice of discarding trajectories 
with large prefactors 23|, |25( does not, by itself, solve this problem for the present system 
since, to achieve convergence, it is necessary to eliminate almost all trajectories that do not 
autoionize. 

To reduce the severity of this difficulty we observe that it is the phase of the pref- 
actor, rather than its modulus, that is primarily responsible for creating the interfer- 
ence in the phase space integration which results in quantization. Thus, as other works 



demonstrate. 
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36 



371 ] it is possible to improve numerical convergence of the autocor- 



relation expression (without adversely affecting the computed energies) by replacing the 
modulus of the prefactor with a less rapidly increasing function of time, while retaining 
the accurate value of the phase. The choice for this replacement is not very critical, but a 
constant value for |-R qp t| is found not to work well here since the rapidly varying integrand 
phases then cause the resulting autocorrelation signal to decay too quickly. In the present 
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work, we have found it useful to replace the computed value of \Rq P t\ with \Rqpt\ 1 ^ 2 - 

Despite this step, we still find it necessary to discard a very small minority of trajectories 
(e.g., about one trajectory out of 10 5 ) having extremely large prefactors (e.g., |-R qP f| > 
10 20 ) which would otherwise destroy the numerical convergence of the calculations. We 
emphasize that the above steps were used only for the calculations of the eZe configuration. 
For calculations of the nonchaotic Zee configuration, the correct HK expressions for the 
prefactors were used and no trajectories were discarded due to large |-R qP i|, 

One final simplification adopted here was to filter out trajectories with energies suffi- 
ciently far from the average energy of the initial state and the target energy eigenvalues. 
Thus, trajectories with positive energies (above the double ionization limit) and sufficiently 
negative energies were eliminated from the phase space integrations. Removal of trajecto- 
ries with low energies was especially important for the eZe configuration since they become 
highly unstable at relatively short times, leading to very large prefactors, and making it 
difficult to achieve convergence for the required time durations. 



IV. RESULTS 



As in the case of the three dimensional helium atom, the eZe and Zee collinear systems 
have bound, resonance, and continuum states. We denote bound and resonance states with 
the symbol N n where N = 1,2,... and n = N, N + 1, . . . are, respectively, approximate 
quantum numbers for the inner electron and outer electrons. For each N, the states N n 
form a Rydberg series converging, as n — > oo, to an ionization threshold with energy E = 
—2/N 2 a.u. The series associated with N = 1 (E < —2.0 a.u.) are bound states while those 
associated with N > 1 are resonances. 

Figure [T] shows the autocorrelation function and spectrum obtained for the eZe config- 
uration using an initial state with qQ = (2.0,2.0), Pq = (0,0) and 7i = 72 = 1-0. This 

I I 

same state was previously used by Simonovic[38] in a quantum mechanical calculation of 
the autocorrelation function for this system. The spectrum obtained from that calculation 
had three visible peaks, corresponding to energies of the ground state li, the first excited 
state I2, and the lowest resonance state 2 2 . The present semiclassical calculations with this 
initial state were performed using 2.3 x 10 7 trajectories having energies between -3.5 and 
a.u. The convergence of the resulting autocorrelation function can be judged by comparing 
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the heavy curve in the left panel of figure [T] with the light curve, which was obtained using 
half the number of trajectories. The harmonic inversion spectrum shown in the right panel 
of figure Q] has peaks for the same three levels observed in the split propagator quantum 
Also shown in the figure are vertical lines denoting quantum mechanical 



calculation. 



energies, computed by the complex rotation method (see Appendix A for details). The 
semiclassical-quantum agreement is seen to be good, even for the ground state 

Figure [5] presents further results for the eZe configuration, obtained using an initial state 
with parameters qQ = (30.0,30.0), Pq = (0,0), and 7i = 72 = 0.1. This state projects onto 
resonances which form the lowest members of the Rydberg series with N = 5 and N = 6. 
Well-converged calculations of the autocorrelation function were obtained using 3.6 x 10 5 
trajectories with energies between -4.0 and a.u. No trajectories were discarded due to 
large prefactors. Agreement with quantum values is seen to be very good. 

Results from these and additional, similar, calculations for singlet and triplet states of 
the eZe configuration are displayed in Table [B The semiclassical values reported are the real 
parts of the harmonic inversion parameters — Ef~ of (1241) . Perhaps due to the reduction in 
the prefactor modulus, the imaginary parts of these parameters are usually about an order 
of magnitude larger than the quantum values, and are not reported here. 

In our treatment of the Zee configuration, lowest energy levels corresponding to n = N 
and n = N + 1 for specific Rydberg series N were calculated by choosing the (q , Po) 
parameters of the initial states to coincide with the location of the central frozen planet orbit 
for this system. Since this orbit passes near the phase space point qo = (2.15,5.76), = 
(0,0) at E — — 1,[3] the initial state for calculation of quantum states at an approximate 
energy E was chosen to have parameters q^ = (2.15/|J5|, 5.76/|i£|), Pq = (0,0). Figure [3] 
shows results obtained in this way with the estimate E = —0.257 a.u. and 71 = 72 = 0.4. 
About 2.3 x 10 6 trajectories with energies between —3.0 and a.u. were used to produce 
the convergent autocorrelation curve shown in this figure. 

Results for these and other Zee states, obtained in a similar way, are reported and 
compared to quantum energies in Table [III Since, for this configuration, our semiclassical 
calculations yielded numerically identical energy levels for symmetric and antisymmetric 
states with n > N, only one set of energies is reported for each set of quantum numbers. 
It appears likely that the singlet-triplet splitting in these cases is due to tunnelling through 
the electron-electron repulsion barrier. It is known that the HK method is unreliable for 
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treatment of such "deep" tunnelling phenomena 



V. DISCUSSION AND CONCLUSION 

Periodic orbit methods are superior tools for understanding how quantization emerges 
from classical behavior. However IVR techniques provide a way for performing semiclassical 
calculations that may remain practical even when the periodic orbit treatments are unfea- 
sible. The results of this work show that the such methods can be successfully applied to 
systems with more than one electron. 

Such applications are, however, not trivial. The Coulomb singularity, the high degree of 
chaos, and the instability of the classical dynamics characterizing these systems, significantly 
complicate the use of IVR methods. Nevertheless, these difficulties can be largely overcome, 
allowing calculation of converged autocorrelation signals for times that are long enough to 
quantize the systems treated. Among the methods used to achieve this were regularization 
of the motion, reduction of the modulus the HK prefactor, removal of trajectories with 
sufficiently negative energies, and application of the harmonic inversion method to form the 
autocorrelation spectrum. 

The problems mentioned in the Introduction, concerning the classical autoionization of 
many-electron systems, turned out to be somewhat less serious than anticipated. The clas- 
sical escape of electrons does require the calculation of additional trajectories, but the com- 
putational effort involved is partly offset by savings obtained by terminating trajectories 
once they ionize. 

Agreement between the semiclassical results and quantum energies was generally very 
good. The semiclassical energies obtained for almost all bound and resonance states exam- 
ined were within about 1% of the exact values. This accuracy was achieved even for the 
ground state of the eZe configuration. We note that previous semiclassical treatments of 

n q 

this state [6|, |7[ provide the current semiclassical estimates of the ground state energy for the 
three-dimensional atom. 

In practice, even with the use of harmonic inversion, resolution of individual levels in the 
energy spectrum requires calculation of the autocorrelation function for longer times when 
the density of states increases. Since, for each N, this density becomes infinite as n — > oo, 
IVR methods are not suitable for calculating levels with large values of n — N . However, 
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as the energy increases (\E\ — > 0), the scaling relation, (fTTl) . implies that IVR calculations 
can be continued for longer times without encountering significant problems due to chaos, 
autoionization, etc. The results obtained here demonstrate that this allows application of 
the IVR method to determine energies for the first few resonances of each Rydberg series 
N, even when N is moderately large. 

As stated in the Introduction, the present work is motivated by the possibility of perform- 
ing semiclassical calculations for the fully three-dimensional helium atom and other atomic 
and molecular systems. Since the number of trajectories needed in the present calculations 
is already large, attempts to treat such systems would certainly benefit from improvements 
in the efficiency of the IVR technique. One way in which this can be accomplished is by 
exploiting the energy scaling of the classical motion more thoroughly. In the present work, 
separate trajectories were calculated for each physical initial condition sampled in the four- 
dimensional phase space. With an appropriate reformulation, it should be possible to achieve 



savings by running trajectories at a single energy only and scaling the results. 22j] A second 



way to increase the efficiency of the calculations is to replace integration over one of the 



remaining phase space variables with an integration over time along the trajectories [39|, |40|. 
Calculations for other systems show that this approach significantly reduces the number of 
independent phase space points which need to be sampled to obtain converged results. Still 
greater improvements can be often be achieved by applying a family of IVR methods which 



cut the dimensionality of the phase space integration by half. 41] 

Even without these modifications, the present HK treatment may prove to be successful 
for many states of the three dimensional helium atom. Since the motion perpendicular to 
the collinear configurations is known to be stable, it may perhaps be said that the present 
calculations for the eZe configuration already treat the most problematic subspace of the 
three-dimensional system. Thus, application of IVR methods seem to place full- dimensional 
semiclassical calculations of atoms and molecules within reach. 
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APPENDIX: QUANTUM CALCULATIONS 



Energy levels for singlet states of the eZe configuration were previously calculated quan 



turn mechanically by Blumel and Reinhardt[4| using the complex rotation method. 



44 



42 



43|, 



45J In addition, three levels for this system were estimated by Simonovic[38] from a 
computed quantum autocorrelation spectum. To obtain values for comparison with our 
calculations, which also include energies for triplet states of the eZe configuration as well 
as energies for states of the Zee configuration, we carried out quantum complex rotation 
calculations that are described in this Appendix. The energies obtained in this way for the 
singlet eZe states, reported in Table [H agree with those of|j. 

Since resonance widths were not reliably obtained with the semiclassical method investi- 
gated in this paper, no effort was made to determine imaginary parts of resonance energy 
eigenvalues in the present quantum calculations. Additionally, since singlet-triplet splittings 
were not resolvable for the Zee configuration by our semiclassical method, our quantum 
calculations for this configuration were limited to singlet states. 

The basis set used for bound and resonance states of the eZe configuration consisting of 
the functions 

4>nm{ri,r 2 ) = N nm [(an) ra (ar 2 ) m ± (an) m (ar 2 )"] e~ a ^\ n,m>l, (A.l) 

where the + and — signs apply for singlet and triplet states, respectively. N nm is a normal- 
ization constant and a = \a\e~ ie is a complex scaling parameter. For Zee singlet states, the 
basis consisted of the functions 

0n m (ri,r 2 ) = N nm (a r< ) n (a r> - a r< ) m e - ar > n,m > 1, (A.2) 

where r< and r> denote, respectively, the smaller and larger of r± and r 2 . In all cases, 
the required overlap and Hamiltonian matrix elements can be readily expressed in terms of 
elementary functions. 

To orthogonalize the basis set and remove linear dependences, the overlap matrix was 
diagonalized and the Hamiltonian was represented in the basis of the orthogonal eigenvectors 
associated with eigenvalues greater than a cutoff value, typically chosen as 1.0 x 10~ 16 . The 
Hamiltonian matrix was then diagonalized with variable basis sizes and fixed complex a 
until convergence was obtained for the real parts of desired eigenvalues. Typical values for 
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the scaling parameters were 9 = 30 rad and \a\ = 2/N (for eZe states) or \a\ = l/N (for Zee 
states). Fewer than 600 basis functions were usually required to achieve convergence for the 
modest precision required here, although for some high-energy resonances, more than 1400 
functions were used. \a\ and 9 were then varied until the real parts of the desired eigenvalues 
were stationary. With sufficiently large basis sets, this condition was often already achieved 
with the initial choice of these parameters so that further optimization was unnecessary. 
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TABLE I: Binding energies 


-Re(E Nn 


)] for symmetric (sin 


iglet) and antisymmetric (triplet) states 


of the of the collinear helium atom in the eZe configuration. SC: values calculated with the present 


IVR technique; QM: energies obtained 


using the quantum mechanical complex rotation treatment 


described in the Appendix. 








State 




Symmetric 


Antisymmetric 


N n 


SC 


QM 


SC QM 


1 1 


3.2102 


3.2459 




2 


2.2225 


2.2028 


2.2393 2.2254 


2 2 


0.8216 


0.8224 




3 


0.6045 


0.6098 


0.6150 0.6184 


3 3 


0.3621 


0.3662 




4 


0.2883 


0.2890 


0.2906 0.2935 


5 


0.2615 


0.2603 


0.2616 0.2618 


4 4 


0.2050 


0.2061 




5 


0.1690 


0.1695 




6 


0.1528 


0.1532 




7 


0.1437 


0.1441 




5 5 


0.1314 


0.1320 




6 


0.1112 


0.1117 




7 


0.1006 


0.1016 




6 6 


0.09127 


0.09204 




7 


0.07922 


0.07928 
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TABLE II: Binding energies for states of the collinear helium atom in the Zee configuration. 
SC: values calculated with the present IVR technique; QM: energies obtained using the quantum 
mechanical complex rotation treatment described in the Appendix. 



State -Re(£^7v n ) 



N 


n 


SC 


QM 


1 


1 


2.1265 


2.1084 




2 


2.0556 


2.0504 


2 


2 


0.5408 


0.5394 




3 


0.5260 


0.5240 


Q 




Q 






n i a on 
U.z4zU 




4 


0.2368 


0.2360 


4 


4 


0.1370 


0.1368 




5 


0.1342 


0.1339 


5 


5 


0.08788 


0.08783 




6 


0.08620 


0.08618 
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-3 -2 

energy (a.u.) 



FIG. 1: Results obtained for the eZe configuration with the initial state having parameters = 
(2.0, 2.0), Pq = (0, 0) and 71 = 72 = 1.0 a.u. The left panel shows the real part of autocorrelation 
functions calculated using about 2.3 x 10 7 trajectories (heavy curve) and 1.1 x 10 7 trajectories (light 
curve). The right panel shows the energy spectrum obtained by the harmonic inversion method. 
The vertical dashed lines indicate quantum levels. 
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time (a.u.) energy (a.u.) 



FIG. 2: Results obtained for the eZe configuration with the initial state having parameters qQ = 
(30.0,30.0), Pq = (0,0) and 71 = 72 = 0.1. The left panel presents two curves for the real part of 
autocorrelation function that are difficult to resolve on the scale of the figure: a thick curve showing 
results obtained using about 3.6 x 10 5 trajectories and a thin curve showing results obtained with 
about half this number of trajectories. The right panel shows the energy spectrum obtained by 
the harmonic inversion method. The vertical dashed lines indicate quantum levels. 
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1000 1500 
time(a.u.) 



0.25 -0.24 -0.23 
energy(a.u.) 



-0.22 



FIG. 3: Results obtained for the Zee configuration with the initial state having parameters 



q<5 



(8.37,22.41), = (0,0) and 7l = 72 = 0.4. The left panel shows the real part of the 



autocorrelation functions. The inset shows details of the real part and modulus of c(t) obtained 
with about 2.3 x 10 6 trajectories (heavy curve) and 1.2 x 10 5 trajectories (light curve). The right 
panel shows the energy spectrum obtained by the harmonic inversion method. The vertical dashed 
lines indicate quantum levels. 
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